Integrating Random Matrix Theory Predictions with Short-Time Dynamical Effects 

in Chaotic Systems 
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We discuss a modification to Random Matrix Tiieory eigenstate statistics, that systematically 
takes into account the non-universal short-time behavior of chaotic systems. The method avoids 
diagonalization of the Hamiltonian, instead requiring only a knowledge of short-time dynamics for a 
chaotic system or ensemble of similar systems. Standard Random Matrix Theory and semiclassical 
predictions are recovered in the limits of zero Ehrenfest time and infinite Heisenberg time, respec- 
tively. As examples, we discuss wave function autocorrelations and cross-correlations, and show 
that significant improvement in accuracy is obtained for simple chaotic systems where comparison 
can be made with brute-force diagonalization. The accuracy of the method persists even when the 
short-time dynamics of the system or ensemble is known only in a classical approximation. Fur- 
ther improvement in the rate of convergence is obtained when the method is combined with the 
correlation function bootstrapping approach introduced previously. 

PACS numbers: 05.45.Mt, 03.65.Sq 
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I. INTRODUCTION 

The statistical behavior of chaotic wave functions has 
been a key topic of investigation from the early history 
of quantum chaos and wave chaos physics, and its study 
is essential for improved understanding of resonances, 
transport, and long-time dynamics in non-integrable sys- 
tems [1| . The standard method of determining the statis- 
tics of the eigenstates and eigenvalues of any physical 
Hamiltonian is by direct diagonalization. On the other 
hand, the semiclassical ergodic hypothesis suggests 
that the statistics of chaotic eigenstates in the semiclas- 
sical limit are universal: just as a typical trajectory in 
a classically chaotic system uniformly covers the entire 
energy hypersurface in phase space, independently of the 
detailed dynamics, so should a typical eigenstate in the 
corresponding quantum system be fully spread over the 
available Hilbert space. Such behavior follows naturally 
from Random Matrix Theory (RMT) Q , which describes 
a statistical ensemble of Gaussian-random Hamiltonians 
having no preferred basis. Within the RMT approxi- 
mation, eigenstates are random vectors either in the full 
Hilbert space or in the subspace given by energy and 
other conservation laws. For a quantum particle in a 
slowly-varying potential, a typical wave function then be- 
haves locally like a random superposition of plane waves 
of fixed wave number, as discussed by Berry RMT 
wave function statistics have been well confirmed as a 
good zeroth-order approximation for a variety of physi- 
cal systems in the absence of integrability. 

The related "quantum chaos conjecture" 043 which 
states that the spectra of individual classically chaotic 
systems also obey RMT statistics in the semiclassical 
limit, is similarly well established experimentally and nu- 
merically. Recently, significant progress has been made 
in deriving this result analytically, starting form a peri- 
odic orbit representation Q. 

As a universal theory, RMT specifically excludes 



system-specific behavior associated with short-time dy- 
namics, boundary conditions, or interactions. Well- 
recognized deviations from random wave function statis- 
tics are associated with boundary effects [1, [13], finite 
system size 9| , unstable periodic orbits [lli > diffusion [l| , 
and two- body random interactions in many-body sys- 
tems (l^ . IT3| . Similar deviations from RMT spectral 
statistics have also been long recognized IJ], and are 
known to arise from non-universal short-time dynam- 
ics d5]. 

Much progress has been made in understanding such 
deviations in specific situations of physical interest, for 
example chaotic wave function correlations in Husimi 
space associated with classical dynamics [l^ and real- 
istic mesoscopic S-matrices arising from a diffusive ray 
picture of wave propagation [l7j . In particular, semiclas- 
sical methods [l8| have proven successful in quantifying 
the effects on wave functions of boundaries d, and 
periodic orbit scars Tl']. However the limit implied by 
semiclassical approximations may not always be achiev- 
able or relevant in describing actual experiments. For ex- 
ample, an analysis of electron interaction matrix elements 
in ballistic quantum dots shows that even for thousands 
of electrons in the dot, several important statistical quan- 
tities typically exceed random wave predictions by a fac- 
tor of 3; for other quantities the ra,ndom wave model even 
fails to predict the correct sign [ilj (see also ^6]). This 
failure of the random wave picture has important impli- 
cations for interpretation of conductance experiments in 
the Coulomb blockade regime. 

In some situations, e.g., [l9j . brute force diagonaliza- 
tion of the Hamiltonian may be used to obtain correct 
statistics for the stationary or long-time behavior, but 
for very large Hilbert spaces, such as those that arise in 
many-body situations, N x N matrix diagonalization is 
often impractical, due to the 0{N^) scaling. Even where 
it "works" , diagonalization is unlikely to produce any in- 
tuition about the relevant physics, and must be repeated 
for each new Hamiltonian. Perturbation theory addresses 
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this situation, but only when the change in the Hamilto- 
nian is small compared to the mean level spacing. In fact, 
individual eigenstates of a chaotic Hamiltonian are highly 
sensitive to perturbations of the system, particularly in 
multi-particle systems. The statistics of such systems 
are much more robust and remain accurate for moderate 
perturbations. For these reasons, a method of predict- 
ing wave function statistics for a family of Hamiltonians, 
similar in form but containing non-trivial perturbations, 
is highly desirable. 

Here we present a system- and basis-independent 
means of supplementing RMT with short time dynamical 
information, an area in which little has been done sys- 
tematically despite it being well known that RMT ignores 
dynamical effects. Our method, introduced previously 
in abbreviated form 
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eliminates the need for diag- 
onalization of the Hamiltonian, while providing greatly 
improved accuracy over standard RMT and semiclassi- 
cal methods for finite systems with a finite Ehrenfest 
time. Instead of treating RMT and non-universal short- 
time behavior separately, we show that they can be natu- 
rally combined to produce useful quantitative predictions 
about wave function statistics in realistic non-integrable 
systems. 

Below, we apply our techniques to eigenstate intensity 
auto-correlations and two-point correlations, as they are 
the simplest nontrivial moments of the eigenstates. Phys- 
ically, these represent the average return probability for a 
given initial state and the average transport probability 
between a pair of states. Our method is easily extendable 
to other moments of the eigenstate probability distribu- 
tion. Similarly, while we focus on the behavior of eigen- 
states in the position or Gaussian wave-packet (Husimi) 
representation, the method applies equally well to eigen- 
state statistics in the momentum basis (of interest for 
scattering problems) or in any other representation. 

A useful way to think about the relationship between 
real non-integrable dynamical systems and RMT is that 
in a real dynamical system, a generic initial wave packet 
spreads (e.g., exponentially with a rate given by the max- 
imal Lyapunov exponent A if the Hamiltonian is chaotic), 
and eventually covers the entire available Hilbert space 
of dimension N. This occurs after the Ehrenfest time 
Te ~ IniV for a ballistic chaotic system. For a dif- 
fusive dynamics, the analogous time scale is the Thou- 
less time Tth (the time required to spread diffusively over 
the system [22|), which is proportional to the square of 
the sample size. RMT, on the other hand, assumes that 
this randomization takes place instantaneously, and the 
behavior is universal for all times t > 0. Thus, RMT 
corresponds to a non-integrable system with an infinite 
Lyapunov exponent (A — oo), or zero Ehrenfest or Thou- 
less time. 

A Lyapunov exponent of infinity, as implied by RMT, 
represents an explosion. The initial state is instantly de- 
stroyed and randomly redistributed across the entire sys- 
tem after a single step regardless of the the details of the 
system. Clearly this is not physically realistic, but it is 



often a useful approximation for several reasons. Firstly, 
for highly chaotic systems, A is large, and the RMT pre- 
dictions are quite accurate. Secondly the predictions of 
RMT for many quantities are computationally simple. It 
is often preferable to accept the errors of RMT in order 
to avoid the computational difficulty of brute force di- 
agonalization. Thirdly, in the case of spectral statistics, 
RMT fails only on large energy scales (corresponding to 
times shorter than Te or Tth), and it is entirely justified 
to apply RMT inside energy windows of size less than the 
Thouless energy Tc ~ fi/Tth (e.g., in the study of level 
spacing statistics). The situation with wave functions is 
not so simple, as short-time dynamical effects leave an 
imprint at arbitrarily long times, and inevitably lead to 
large deviations from RMT not only for energy-averaged 
quasimodes but also for individual eigenstates [231 . 



II. METHOD 

To enable direct comparison with RMT, let us con- 
sider fully chaotic (ballistic or diffusive) dynamics with- 
out symmetry on an A^-dimensional Hilbert space with 
eigenstates |^„). To avoid ambiguities in the definition 
of l^n), we assume a non-degenerate spectrum. Typical 
quantities of interest, then, are functions of the ampli- 
tudes (a|^„) for any physically-motivated basis state |a), 
which may be a position or momentum state, a Slater 
determinant for many-body systems, or more generally 
an eigenstate of some zeroth-order Hamiltonian. Given 
the normalization condition X^^iLi = 1' ^^'^ sim- 

plest and first non-trivial moment of these amplitudes 
is given by the local inverse participation ratio (IPR), 
which measures the degree of localization at \a): 
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(1) 



varying from P°° = 1 in the case of perfect delocalization 
to P"" = N for perfect localization. For two arbitrary 
states we similarly define 
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P''' = Nj2\{a\CMb\U\\ 



(2) 
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a measure of intensity-intensity correlation. These quan- 
tities are directly related to the time-averaged return 
probability and transport probability, 

P^^=N\im ^ r dt\{a\am' (3) 



and 



P'''=N\im ±, r dt\{a\bm\ (4) 

T — >cxD Zl J _x 

as may be easily verified by inserting two complete sets 
of states J2n in Eq. © or dH). 
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In the more general case of a time-periodic Hamilto- 
nian H (t + Tg) = H(t), to be considered in the numeri- 
cal examples below, eigenstates |^„) will be the Floquet 
eigenstates, |^n(?o)) — e~'"^""^°/''|^„), and the time inte- 
grals in Eqs. dSj and (|4]) are replaced by discrete sums 
with time step Tq. 

Obviously, higher-order moments of the eigenstate in- 
tensities and in general the entire joint distri- 
bution of the intensities, may be considered (e.g., 
in similar fashion. We may also relax the requirement 
that only pure states such as \a){a\ act as probes, and in- 
stead measure the structure of chaotic eigenstates using 
any desired self-adjoint operator a [1^. Operator probes 
(of phase space size greater than or smaller than h) will, 
for example, be particularly helpful in the study of hi- 
erarchical eigenstates in a mixed chaotic-regular phase 
space [11]. Again, without loss of generality we may 
adopt the normalization Tra ~ 1. Eqs. and (|3]) be- 
come 

N 



= ^ lim ^ I dtTram, (5) 

with the autocorrelation P"" as an obvious special case. 

In the limit N — >■ oo, where the Heisenberg time at 
which eigenstates are resolved becomes long compared 
to all other relevant time scales in the problem, averages 
of the form (H}-® may be obtained using short-time dy- 
namics. [In the following, we refer to this limit as the 
"semiclassical" limit, whether or not the system Hamil- 
tonian has a classical analogue.] Specifically, using the 
correct normalization for discrete-time dynamics we have 
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where have defined the time-domain quantity 

P'^\t)^\{am)\' + {<^W))m\h), (7) 

and where \a) and \b) may be any two states (identical, 
orthogonal, or in general overlapping). The \a) = \b) case 
is the special case of the IPR and we have 



P"''(i) = 2|(a]a(i))p; 
thus Eq. (IH) becomes 



P- «2 ^ \{a\am 



(8) 
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time T must be long compared to the ballistic or diffusive 
Thouless time [13] (so as to include all the non-universal 
dynamics), and short compared to the Heisenberg time, 
which scales with the effective Hilbert space dimension 
N. 

Importantly, no distinction is made in Eq. ([B]) between 
non-universal short-time revivals that indicate deviations 
from RMT in the eigenstate statistics and the 0(1/A^) 
short-time revivals that are present already in the context 
of RMT. Effectively, these random revivals are double- 
counted as both short-time and long-time returns. As 
a result, Eq. ^ systematically overestimates corrections 
to RMT, and violates probability conservation P"^ = 
1 given a complete basis \b) for any r > 0, with the 
violations growing linearly as t/N. (In particular, the 
violation necessarily disappears in the limit iV — ^ oo, 
where a clean separation of scales exists between short- 
time and long-time behavior.) 

We now notice that all problematic aspects of Eq. © 
for finite system size N can be eliminated by introducing 
a T- and (a]&)-dependent prefactor: 
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(10) 



where in particular C^'^^(t) 



N/At converges to the 
oo (similarly to Eq. (jH). To fix 

C)^''''(t) for general t, we apply Eq. ([TOl) to the case of 



exact answer as r 

(q|6) 



RMT dynamics and obtain 
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Eq. ([TT|). and its generalizations, constitute the key result 
of this paper as they demonstrate the method of combin- 
ing short time dynamics with RMT to produce accurate 
predictions for eigenstate statistics without diagonaliza- 
tion. It is natural to extend this approach to higher-order 
moments, for example. 







(12) 



^/-r dt P^Mri^) J 

The results generalize also to the statistics of operator 
expectation values, e.g., P"^ and (P"^)^, where P"^ is 
defined by Eq. ([5]). For example. 



_ r dtP°"^{t) 



RMT 



(13) 



In Eq. ([SJ and in the following, indicates an en- 
semble average. If desired, the ensemble of Hamiltonians 
may be selected so that all realizations possess the same 
short-time dynamics P"'^{t), in which case the average 
on the right hand side of ([5]) is superfluous. The cutoff 



The predictive power of Eqs. ([TTl), and ^ will 

be demonstrated in Sec. IIVI below. 

Thus, stationary eigenstate properties of a quantum 
chaotic system or ensemble of systems may be fully de- 
scribed by a combination of short-time dynamics for that 



system or ensemble, with exact results from RMT, with- 
out any need for matrix diagonalization. Reassuringly, 
Eq. (fTTj) and its generalizations yield exact results in 
three limits of interest. In the RMT limit, i.e., for a sys- 
tem with infinite Lyapunov exponent, we have P°''^{t) — 
-Prmt(^) and thus Eq. ^ predicts P"'' = P^mt Sec- 
ondly, in the semiclassical limit N/t — )■ oo, we find 

-Prmt ~ S-r^RMT(0 (i^ discrete-time case), and 
from Eq. we recover Eq. Finally, in the limit 
where an infinite amount of dynamical data is available 
as input, r — >■ cx), the prediction of Eq. (1111) becomes 
exact due to the identity (HJ. 

More importantly, as we will see in Sec. IIV| Eq. (fTTj) 
and its extensions provide reliable approximations to ex- 
act diagonalization in situations quite far from any of 
these three limits, i.e., for finite-sized weakly chaotic sys- 
tems, and where the only input is short-time dynamics 
on the scale of a few Lyapunov times. 



III. EXPLICIT EXPRESSIONS: BROKEN TIME 
REVERSAL SYMMETRY 



appropriate to that ensemble, e.g 
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V (e«(£„' (16) 
^ V /RMT 



(l(«l?n)PI(&l?n'>P + {a\in){Uh){^n'\a)m,) 



RMT 



For the specific case of discrete-time Floquet dynamics, 
described by the CUE ensemble, which will be relevant 
for the numerical examples in Sec. IIV[ we have 
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for t = 
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N+l 



P^lATit) - (1 + \{a\b)\') X <! i±^^ for 1 < |t| < iV , 

for |t| > N 

(17) 

and analogous expressions for self-adjoint operators are 
obtained by spectral decomposition, as in ((T5]). Results 
for other RMT ensembles (Gaussian or circular) are ob- 
tained by evaluating Eq. ([T6l) with the spectral form fac- 
tor and eigenstate statistics appropriate to that ensem- 
ble. 

For the discrete-time case, Eq. pO)) takes the form 



The short-time overlaps P°^(i) needed as input to 
Eq. (jlip may sometimes be known analytically, as in the 
case of periodic orbit scars, while in more general sit- 
uations the short-time dynamics for a given system of 
interest is obtainable numerically, to any desired time 
scale T. The RMT factors in PT|) and its generalizations 
will depend on the specific random matrix ensemble of 
interest (e.g., orthogonal, unitary, or symplectic; Gaus- 
sian or circular). For a given random matrix ensemble, 
these factors may be calculated numerically or treated 
entirely analytically. 

For example, in the absence of time reversal symmetry 
(i.e., for the GUE or CUE ensembles), we have standard 
results (for arbitrary \a) and \b)), 
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^RMT — ^ 
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(14) 



while for a general self-adjoint operator a we apply spec- 
tral decomposition and obtain 



P<^^ « Cj;l^'(r)^P°''(i). 



(18) 



Now, substitution of Eqs. ([Ml) and dTT]) into Eq. ^ 
yields the explicit closed-form result 



pab r^^pab(l-^ 



for < T < iV 



f N 

X J N+1+2t+(t-^+t)/N 

X^TTr forr>7V 



(19) 

Here we note that all dependence on the initial-state over- 
lap (a 1 6) cancels in the final expression (fTO]) . so the final 
result has the same form for identical, overlapping, or or- 
thogonal states |a), \b). As TV — > oo for fixed r, Eq. (flQ)) 
reduces to the semiclassical expression ([6]), while in the 
limit r ^ A'^ of complete dynamical information, Eq. (|19p 
reduces to the exact answer, 



N ^ 

pab ^ lij^ _\ pabU\ 
r-i-oo 4t'^ ' 



(20) 
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N + l 



(15) 



where Ai are the eigenvalues of a (Yin Ai = 1). Similar 
expressions exist in the presence of time reversal symme- 
try, where one must be careful to distinguish the situ- 
ations of time-reversal invariant and non-invariant basis 
states (27| . 

Similarly, RMT dynamical overlaps may be expressed 
exactly for any random matrix ensemble using RMT 
eigenstate statistics and the RMT spectral form factor 



the discrete-time form of the identity (|4]) . 



IV. MODEL AND NUMERICAL RESULTS 
A. Model 

We now discuss a few illustrative examples, using as 
our model the paradigmatic example of a quantized peri- 
odically kicked Hamiltonian [l^l on the compact toroidal 
phase space (g,p) G [-1/2, 1/2) x [-1/2, 1/2). To enable 
direct comparison with RMT, we consider fully chaotic 
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(ballistic or diffusive) dynamics without symmetry. The 
time-periodic Hamiltonian is given by 



H{q,p,t)=Tip) + V{q) ^{t - nT^) 



(21) 



where without loss of generality we will set the period Tq 
to unity. Evolution for one period defines a discrete-time 
map 



p{t + l)=p{t)~V'{q{t)) 

q{t + l) = q{t)+T'{p(t + l)), 



(22) 
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where we have chosen to begin each period immediately 
before the kick. The kinetic term T{p) and potential kick 
term V(q) are selected to produce a perturbed Arnold's 
cat map [29|, 

T{p) = + ^ cos(27rp) -I- t{p) 
V{q)^~'^q^-^COs{27Tq)+v{q), (23) 

where the original (linear) Arnold's cat map is obtained 
for K = t{p) — v{q) = 0. In the following, we take 
\K\ < m to guarantee full chaos (for t{p) — v{q) — 0). 

It is easy to see that the map ((22)) has a fixed point 
(i.e., a periodic orbit of period 1) at g = p = 0. Since the 
dynamics is fully chaotic, this fixed point is necessarily 
unstable, and in fact a straightforward calculation shows 
that its instability exponent (local Lyapunov exponent) 
is Aq = cosh~^ (l -f (m — K)^/2j, which reduces to Xq ~ 
m — K for m — K <^ 1. Thus, the perturbation parameter 
K, in addition to breaking the linearity of the original 
cat map, makes it easy to control the chaoticity of the 
dynamics. 

Finally, the perturbations t{p) and v{q) in Eq. ([23]) 
serve to break the parity and time reversal symmetries of 
the Hamiltonian. We choose these functions to be ran- 
dom within a small region near the edges of the phase 
space {\p\ > 1/2 — A and |g| > 1/2 — A) and zero else- 
where. In the following, we set A = 0.1, but the results 
have no significant dependence on A. Conveniently, the 
need for symmetry breaking simultaneously allows us eas- 
ily to construct ensembles of similar Hamiltonians (hav- 
ing the same local Lyapunov exponent Aq), by varying 
t(p) and v(q) in Eq. (j23p while keeping m and K fixed. 

Quantization of the discrete-time dynamics (|22l) on 
a Hilbert space of dimension N is given in the usual 
way i29,] by the one-step unitary evolution operator 

U = g-^f{p)/h^~^V{q)/h ^ ^24) 

where h — 1/2ttN, q and p both take on discrete values 
j/N {-N/2 < j < N/2), and the eigenstates of U are 
the Floquet eigenstates |^„). 



FIG. 1: (Color online) The inverse participation ratio P"" 
for a Gaussian distribution centered on a short periodic orbit 
with local Lyapunov exponent Aq — 0.5 is computed by direct 
diagonalization, solid lines, and compared with the short-time 
dynamical prediction given by Eq. squares. Here the 

system size is A'^ = 32 and the Gaussian distribution has size 
s = 0.5 (Left panel) or s = 0.25 (Right panel). Convergence 
of the dynamical prediction to the exact result is observed 
when the dynamical calculation includes information about 
times r up to 2 or 3 in units of the local Lyapunov expo- 
nent Ao. The RMT value Prmt = (1 + s~^)N/{N + 1) and 
the semiclassical result Pg^ = (H- s~^) X]t^-oo sech(Aot) are 
shown for comparison. All quantities appearing here and in 
subsequent figures are dimensionless. 



B. Inverse Participation Ratio 

We begin by considering the inverse participation ra- 
tio P"", where a is the Weyl transform of a Gaussian 
distribution p{q,p) ^ e '"i ^ '"p centered on the pe- 
riodic orbit g = p = 0. We define s = UqUp/h as the size 
of the Gaussian distribution. In the special case s — 1, 
d is a projection onto a minimum uncertainty Gaussian 
wave packet, while more generally a represents a mixed 
initial state, corresponding to a classical phase space area 
smaller than or greater than h. 

Results are shown in Fig. [1] where the dynamical pre- 
diction of Eq. (ITT]) for several values of the cutoff time r 
is compared with exact values obtained by brute-force di- 
agonalization, and also compared with the semiclassical 
and RMT approximations. We note that the dynamical 
prediction is identical to the RMT limit for r = 0, as 
it must be since the ratio of the integrals in Eq. (jlip is 
1 at T = and the method thus reduces to the RMT 
prediction when no dynamical information is available. 
We also notice that the prediction for the inverse partic- 
ipation ratio quickly converges to the exact answer after 
only 2 or 3 Lyapunov times. 

Fig. [21 illustrates the relationship between the exact 
value of the inverse participation ratio, our dynamical 
prediction, and the limiting RMT and semiclassical ap- 
proximations, as the system size TV is varied. The RMT 
approximation fails entirely for all system sizes except 
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FIG. 2: (Color online) The inverse participation ratio for a 
pure Gaussian wave packet (s = 1) is computed exactly and 
compared with the dynamical prediction of Eq. using 
dynamical information up to times r = 2Xq^ and 4A(^^, where 
Ao = 0.25 is the local Lyapunov exponent. Results are shown 
for various values of the system size A'^. The semiclassical 
and RMT limits are also shown for comparison (see Fig. [T] 
caption) . 



the trivial case N = 1 (not shown). The semiclassical 
limit is obtained as N ^ oo, but we note very significant 
deviations from this limit even when N takes values of 
100 or greater. The strong deviations from both RMT 
and semiclassical predictions are well reproduced in our 
dynamical calculation for all system sizes N, even using 
a short calculation time Aqt = 2. 




FIG. 3: (Color online) The inverse participation ratio in the 
position basis for Ao — 0.0625 is computed exactly and com- 
pared with the dynamical prediction of Eq. Results are 
shown for several system sizes A'^ with the short time dynam- 
ics carried out to various times r. 



cal Lyapunov exponent is chosen to have a small value, 
Ao — 0.0625, so that we are very far from the RMT 
limit. Furthermore, examples of relatively small N are 
included, where the system is very far from the semi- 
classical limit as well. As seen in Fig. [3l excellent con- 
vergence to the exact results obtained by diagonaliza- 
tion is attained in only a fraction of a Lyapunov time, 
r w 0.5A,7^ 



C. Classical Dynamics Approximation 

In all the above examples, we have considered situ- 
ations where the quantum dynamics up to some short 
time r are assumed to be known exactly, and can be 
used to calculate properties of the eigenstates with arbi- 
trarily good accuracy as t — oo. In practice, however, 
it is commonly the case that even the short-time behav- 
ior is known only in some approximation. For example, 
we may not have access to the dynamics under the full 
Hamiltonian H but only to the dynamics under some 
simpler Hamiltonian Hq. Alternatively, we may be in- 
terested in the eigenstate statistics for an ensemble of 
Hamiltonians Hq + SH, where the short-time dynamics 
is dominated, but not completely determined, by Hq- Fi- 
nally, in many physical situations of interest, the "true" 
system Hamiltonian is not known at all, e.g., in meso- 
scopic experiments where the details of individual eigen- 
state structure change completely from one identically 
prepared sample to the next. In all these situations, we 
would like to know the robustness of Eq. (fTTj) and its gen- 
eralizations, i.e., how well can we predict the statistical 
properties of eigenstates given only approximate knowl- 
edge of short-time dynamics? 

In Fig. m we consider a situation similar to that 
presented in Fig. [2l but where the short-time dynam- 
ical overlaps P°'°'{t) for < t < t are calculated ei- 
ther exactly, as in Eq. ([5]), or using classical evolution, 
-PclasW = 2jjdqdpp{q,p,t)p{q,p), where p{q,p) is the 
phase space density corresponding to the wave packet \a). 
The predictive power of the classical short-time dynamics 
is almost as good as that of the exact quantum short-time 
dynamics for all system sizes N considered, once correct 
normalization by the RMT factors is implemented, as 
prescribed by Eq. pT|) . In each short-time cut- 

off r equal to twice the local Lyapunov time is used for 
both quantum and classical dynamics; obviously further 
improvement in the predictive power would be obtained 
by using larger r. Once again, both the semiclassical and 
RMT approximations work very poorly for most system 
sizes of interest. 



The applicability of the method is in no way limited 
to Gaussian test states. In Fig. [3] we compare the pre- 
dicted and actual values of the inverse participation ratio 
in position space, J2q ^''j where the sum is over all 
N position states \q). The calculation is performed for 
various system sizes N and short-time cutoffs r. The lo- 



D. Interaction Matrix Element Variance 

As another example, we consider wave function inten- 
sity correlations P"'' for position states |a), \b). Since 
1^ 6=1 = 1 is given identically by wave function 
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FIG. 4: (Color online) The inverse participation ratio for a 
pure Gaussian wave packet (s = 1) is computed exactly and 
compared with the dynamical prediction of Eq. using 
either exact quantum dynamical information P°"(t) or the 
classical evolution -Pc'ias('') uP to times r = 2\q^ . Here Ao ~ 
0.125 is the local Lyapunov exponent, and results are shown 
for various values of the system size A^. The semiclassical 
and RMT limits are also shown for comparison (see Fig. [T] 
caption). 



normalization, we focus on the first interesting moment, 
the variance 



(25) 



,.b=l 



is a sirtiple measure of non-uniformity in infinite-time 
transport [301 , and ranges from Vl^ = for perfect ergod- 
icity to = A'^ — 1 for perfect localization. 

We note also that we can interchange the roles of evo- 
lution eigenstates and basis states, |^„) O |a). W may 
then be equivalently written as the variance of the inter- 
action matrix elements 



N 



(26) 



o=l 



which describe the intensity overlap of eigenstates |^„) 
and in the chosen basis. Note that the summation 
here is over the basis states |a), and not over the eigen- 
states \^) as in Eq. The variance W is then given by 
a summation over all pairs |^„) and |^n'), 



W 



1 

iV2 



N 



(27) 



The statistics of such interaction matrix elements in 
chaotic systems frequently appear in applications ranging 
from quantum dot conductance in the Coulomb blockade 
regime [l9| to control ling directional emission properties 
in microcavity lasers [3l|. 

To predict the variance W using a combination of RMT 
and short time dynamics, we may simply apply Eq. (jl2p 
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FIG. 5: (Color online) The interaction matrix element vari- 
ance W is computed exactly (Eq. (I25p or (|27|) ') and compared 
with the short time prediction (I12|) . with Xqt — 1. Here 
Ao — 0.125. The RMT result and the semiclassical limit 
Wsc ~ l/N are also shown for comparison. At finite sys- 
tem sizes, both RMT and the semiclassical approximation 
fail dramatically. 



to predict (P"'')2 for each pair of initial and final posi- 
tion states \a) and and sum the results. However, 
we now take the opportunity to note that the origi- 
nal intensity correlators P"-'' predicted by Eq. (fTTj) are 
not guaranteed to satisfy the normalization condition 
'W^^b^i^'^^ ~ which holds for the exact correla- 
tors. Exact normalization is only guaranteed in Eq. (jlip 
when r — oo. This in general creates a small overall 
normalization error, which may affect the prediction of 
the higher-order moments from the short time dynamics. 
However, in applying Eq. (fT2]) to predict W, we already 
need knowledge of the short-time dynamics, P'^^{t), for 
every pair of initial and final states |a), \b). Therefore 
with little added computational effort we may correct for 
the normalization effect, further improving the conver- 
gence with r, simply by rescaling by the norm: 



I 1 Y^^f pa'6'^2 



(28) 



Fig. [5] shows that the semiclassical and RMT predic- 
tions are very similar for the system we consider here, 
and both deviate significantly from the exact results for 
finite N . The method presented here, including rescal- 
ing, converges toward the exact answer very quickly, on 
the order of the Lyapunov time, even where the RMT 
and semiclassical predictions are off by a factor of 2 or 
3. The accuracy can be improved further by utilizing the 
dynamics for longer times r. 



E. Bootstrapping 

Our method scales easily to arbitrarily large iV, but 
for extremely large system sizes it may become time con- 
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suming to calculate even a few time steps of the evolu- 
tion. (Note that brute-force diagonalization for a sys- 
tem of such size is likely to be intractable.) To reduce 
the amount of time-domain information needed as input, 
one may utilize the correlation function bootstrapping 
approach 24], which allows approximate evolution dy- 
namics to be extracted starting from as little as one uni- 
tary time step or bounce. Here we briefly review this 
approach, for the simplest case of a single basis state \a) 
and discrete-time dynamics; a more general and detailed 
discussion is found in Ref. 



24 



Let = {a\a{t)) be the return amplitude after 

time t. It is advantageous to partition this total return 
amplitude at a given time t into two parts: the amplitude 
^ncw(0 that is returning to the subspace spanned by \a) 
for the first time at time t and amplitude that has already 
returned to this subspace at times intermediate between 
and t. Att= 1, we have of course = ^""(1). 

Continuing, we have 



A"" (3) 



Aaa 
new 



aa 
new 



A aa 
new 



A aa 
new 



(1) 
(If 
(1)3 

m 



A 

new 
aa 
new 



(2) 
{1)A 
{1) + A 



aa 
new 
aa 
new 



(2) 
(3), 



(29) 



Note that in the expression for A"" (3), the term 
^ncw(l)^now(2) IS associatcd with the portion of the 
wave packet that leaves the initial subspace, returns after 
the second step and remains after the third step, while 
the term A'^'^^{2)A'^'^^{1) is associated with amplitude 
that remains in the initial space after the first time step, 
then leaves and returns again after two additional time 
steps. These two terms are equal in our case of a single 
wave packet, but would be distinct in general. 

For general t we can express the total return amplitude 
as follows, 



E E n^new(i.), 

j=l iiH htj=t fe=l 



A^^it) 
or equivalently 

t 

t' = l 

Inverting Eq. ([31]) . we can write A'^^{t) as 
AZ^{t) - A'^'^it) - ^""(^ - ^')^Sew(i') 



(30) 



(31) 



(32) 



Now assuming the exact returns A°'°'{t) are known as 
before only up to time t, Eq. p2p yields the new retuns 
^now(0 up to time T. We now make the approximation 
that all new returns after time t vanish, v4°°^(i) = for 
t > T, and use Eq. (|5D|) or PT|) to calculate the total 
returns at arbitrary times. Obviously Eq. pop or (j3ip 
produces the exact total returns by construction when 



t < T, while for i > r it is only an approximation. For 
example, in the simplest case r = 1 where only the return 
after one step is available as input, we have A"°(l) = 
A°°(l) (exact), A'"'{2) « (^'^'^(l))^, ^'^'^(3) « (^'^'^(l))^, 
etc. 

This bootstrapping process can be extended to an ar- 
bitrary time f > T, or if desired even to f = oo. Note 
that the more time steps are used as input for the boot- 
strapping, i.e., the larger the base time r is, the more 
accurate the result. 

We have found that the base time r needed in practice 
is often quite small [24], making the correlation boot- 
strapping approach, in combination with our dynamical 
method, computationally viable even for large systems. 
Fig. [S] shows the results for the inverse participation ra- 
tio of a Gaussian wave packet, as in Fig. [21 but for the 
case Aq — 0.5, and using only two time steps of infor- 
mation as input, r = 2. Here the squares represent a 
calculation performed directly using return probabilities 
up to r = 2 in Eq. ([TT|) . while the crosses are calculated 
by bootstrapping these same known recurrences up to a 
longer time f — 2r, and then inserting return probabili- 
ties up to time f into Eq. (jlip . Clearly, dynamical infor- 
mation on the scale of just one Lyapunov time is in gen- 
eral insufficient to accurately describe eigenstate statis- 
tics, but adding bootstrapping significantly improves the 
accuracy, at least for moderate to large system sizes. 
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FIG. 6: (Color online) The ensemble averaged IPR for a 
minimum uncertainty wave packet, as in Fig. [21 but for an 
ensemble of systems with local Lyapunov exponent Ao = 0.5. 
Predictions with and without bootstrapping are compared to 
the exact IPR, to the RMT approximation, and to the semi- 
classical limit. Squares: exact returns up to times t — 2 are 
used as input to predict the IPR using Eq. (Ilip . Crosses: ex- 
act returns up to time r = 2 are bootstrapped (Eq. pop ) to 
extrapolate the dynamics to time t — 4, and this dynamical 
information is then used as input to predict the IPR using 
Eq. dni)- 
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V. CONCLUSION 

We have developed a method that improves on RMT 
eigenstate statistics for non-integrable systems by sys- 
tematically incorporating short-time dynamics. It may 
be applied to general quantum chaotic systems, which 
may or may not have a classical limit. For systems in 
which a classical limit exists, the approach correctly in- 
cludes the effect of non-universal short periodic orbits, 
while simultaneously correcting for the fact that semi- 
classical approximations break down at low energies or 
small system sizes. Thus, we are able to obtain quanti- 
tative predictions both in the semiclassical regime, and 
also in the regime where the semiclassical approximation 
fails. 

The method is conceptually appealing, in that it is 
completely independent of the system being studied 
and it is also independent of the basis of interest. In 
other words, while the method necessarily uses the non- 
universal time-domain behavior of the system as input 
to improve on RMT predictions, it is in no way depen- 



dent on special properties of the system for success. Our 
method is computationally simpler than brute-force di- 
agonalization, and significantly more accurate than RMT 
or the semiclassical limit for realistic systems. 

The rate of convergence can be improved even more 
using the correlation function bootstrapping procedure. 
When exact short-time dynamical information is not 
available, approximate short-time dynamics, e.g., in a 
classical approximation may be successfully utilized. The 
approach can be extended to consider additional eigen- 
state statistical measures beyond those studied here, as 
well as to include symmetry effects (including time re- 
versal symmetry), mixed phase space j32| . and resonance 
wave function statistics in open systems. 
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